#! /USR/BIn/env python
# -*- coding: utf-8 -*-
# vim:fenc=utf-8
#
# Copyright © 2018 crane <crane@crane-pc>
#
# Distributed under terms of the MIT license.

from operator import mul, add, sub
import math
from linear_tools import scale_vector


# def norm(v1):
def length(v1):
    '''
        output: length of v1(向量的模norm/几何长度length)
        可能重命名为length更好些
    '''
    square = lambda x : x*x
    square_sum = sum( map(square, v1) )
    return math.sqrt(square_sum)

def normalize(v1):
    '''归一化模长度到一'''
    unit = scale_vector(v1, 1/length(v1))
    return list(unit)

def make_vector(p1, p2):
    ''' 两个点组成的向量
        point p2 - p1
    '''
    return list(map(sub, p2, p1))

def proj_vector(from_v, to_v):
    '''
        from_v 投影到to_v上, 計算投影後的向量
    '''
    x = dot_product(from_v, to_v) / dot_product(to_v, to_v)
    return scale_vector(to_v, x)

# =================== dot product  =====================
def dot_product(v1, v2):
    '''
    dot_product/inner_product 只支持两个向量, 输出一个标量scalar.
    对应numpy.dot
    '''
    assert len(v1) == len(v2), (v1, v2)
    product_vector = map(mul, v1, v2)
    return sum(product_vector)

def is_perpendicular(v1, v2):
    ''' 两个向量在几何上是否垂直
        或者是否正交: orthogonal(是否互相影响)
    '''
    return 0 == dot_product(v1, v2)

def project_len(vfrom, vto):
    ''' vfrom 在vto上投影的长度 '''
    return dot_product(vfrom, vto) / length(vto)

# =================== cross product  =====================
def cross_product(v1, v2):
    '''
    cross_product, 输出一个向量v3, 垂直于(v1,v2)组成的平面(法向量norm_vector)
    |v3| = |v1 . v2|,  v3的模等于v1,v2组成的平行四边形面积
    目前只支持3个分量的向量.

    outer product, 对应numpy.cross

    这里不需要使用右手定则, 判断生成的向量方向. 右手定则是给人看的.
    这里直接可以计算出准确的向量值, 信息中自然包含方向信息!

    其实是determinant的求解方式
    '''
    assert len(v1) == len(v2)
    a1, a2, a3 = v1
    b1, b2, b3 = v2

    i =   a2 * b3 - a3 * b2
    j = -(a1 * b3 - a3 * b1)        # NOTE: 一定不要忘记取负号
    k =   a1 * b2 - a2 * b1
    return type(v1)([i, j, k])

def area_parallelogram(v1, v2):
    '''
        output: area of parallelogram maked by (v1, v2)
        需要用到叉乘: cross_product
    '''
    return length(cross_product(v1, v2))

def area_triangle(v1, v2):
    '''
        output: area of triangle maked by (v1, v2)
        需要用到叉乘: cross_product
    '''
    return area_parallelogram(v1, v2) / 2

def norm_vector_of_v(v1, v2):
    ''' 求v1向量  v2向量组成平面的法向量 '''
    return cross_product(v1, v2)

def norm_vector_of_p(p1, p2, p3):
    ''' 由point1, point2, point3组成的平面的法向量 '''
    v1 = make_vector(p1, p2)
    v2 = make_vector(p1, p3)
    norm_vector = norm_vector_of_v(v1, v2)

    # distance_to_origin = dot_product(norm_vector, p1) / length(norm_vector)  # 平面上任意一点和法向量的点积(等于平面到原点距离)
    # d2o2 = dot_product(norm_vector, p2) / length(norm_vector)
    # d2o3 = dot_product(norm_vector, p3) / length(norm_vector)
    distance_to_origin = project_len(p1, norm_vector)
    d2o2 = project_len(p2, norm_vector)
    assert distance_to_origin == d2o2
    # print('origin ----> plane distance %s' % distance_to_origin)
    return norm_vector

def distance_to_plane(point, p1, p2, p3):
    ''' point 到(p1, p2, p3)三点平面距离 '''
    norm_vector = norm_vector_of_p(p1, p2, p3)
    connect_vector = make_vector(point, p1) # point 到平面某点的向量(任意选择一点即可)
    # return dot_product(norm_vector, connect_vector) / length(norm_vector)
    return project_len(connect_vector, norm_vector)

# =================== triple_product =====================
def triple_product(v1, v2, v3):
    ''' mixed product (混合积)
        box product(用来求体积, 像一个box)
    '''

    # method1: dot * cross
    cross = cross_product(v2, v3)
    value = dot_product(v1, cross)
    return math.fabs(value)

    # method2: 用求解行列式的方式来算!
    mat = [v1, v2, v3]
    c0 = v1[0] * (v2[1] * v3[2] - v2[2] * v3[1])
    c2 = v1[1] * (v2[0] * v3[2] - v2[2] * v3[0])
    c3 = v1[2] * (v2[0] * v3[1] - v2[1] * v3[0])
    return math.fabs(c0 - c2 + c3)

def volume(v1, v2, v3):
    ''' (v1, v2, v3)组成的平行六面体的体积, 和(v1,v2,v3)顺序无关.
        TODO :这里也可以用(v1, v2, v3): determinant行列式来作为体积
    '''
    return triple_product(v1, v2, v3)

# =================== degree/ arc =====================
def cos(v1, v2):
    cos_theta = dot_product(v1, v2) / (length(v1) * length(v2))
    return cos_theta

def degree(v1, v2):
    cos_theta = cos(v1, v2)
    return math.acos(cos_theta) # arc cos

def determinant(vectors):
    ''' 计算行列式的值 '''
    pass


# =================== main test =====================
def main():
    print("start main")

    print(dot_product([1,0,0], [1,2,4]))

    print( is_perpendicular([1,1,0], [-1, 1, 0]))   # 直角

    print("cross : %s" % cross_product([1,0,0], [0,1,0]))
    print("cross : %s" % cross_product([1,1,0], [-1,1,0]))
    print(cross_product([2,0,0], [0,2,0]))
    print(length([3,3,0]))

    print("area : %s" % area_parallelogram([1, 1, 0], [-1, 1, 0]))    # 2
    print("area : %s" % area_parallelogram([1, 0, 0], [0, 1, 0]))     # 1

    print("triangle area : %s" % area_triangle([1, 0, 0], [0, 1, 0]))     # 1

    print("volume: %s" % volume([1, 1, 0], [-1, 1, 0], [0, 0, 2]))
    print("volume: %s" % volume([0, 0, 2], [-1, 1, 0], [1, 1, 0]))
    print("volume: %s" % volume([-1, 1, 0], [1, 1, 0], [0, 0, 2]))

    print('degree %s' % degree([-1, 1, 0], [1, -1, 0]))

    print('normal_vector %s' % norm_vector_of_p(
        [1, 2, 3],
        [2, 3, 5],
        [1, 5, 4]
    ))

    print('point-plane distance %s ' % distance_to_plane(
        [0, 0, 0],
        [1, 2, 3],
        [2, 3, 5],
        [1, 5, 4]
    ) )

    print('normalize [2,2,2] is %s' % normalize(
        # [1, 2, 3],
        # [1, 1, 1],
        # [0, 1, 0],
        [2, 2, 2],
    ) )


if __name__ == "__main__":
    main()
